SoSe2021
Folienübersicht
(Nicht-parametrische) Mehrstichprobentests
Mittelwerte von ≥ 3 Stichproben vergleichen
Das Problem der multiplen Vergleiche
Beispiel
- Auswirkung der Meerestiefe auf das bakterielle Wachstum
- Probenahme in 5m (Kontrolle), 10m, 20m, 40m, 80m, 160m
→ Warum nicht jede Tiefe mit der Kontrolle mittels t-Test vergleichen?
‘Familywise Error Rate’ (FWER)
- Sehr ineffizient!!! 5 Tests sind notwendig
- Je öfter ein Test mehrfach auf eine Hypothese angewendet wird, desto höher ist das Risiko eines stat. Fehlers 1. Art (\(\alpha\); Ablehnen einer korrekten H0).
- Die FWER ist dabei definiert durch die Wahrscheinlichkeit, mindestens eine der Testhypothesen fälschlicherweise abzulehnen:
- \(\alpha_{FWER} = 1-(1-\alpha)^c \Rightarrow 1-(1-0.05)^5 = 0.23\) (c = Anzahl der Vergleiche)
- => Der gemeinsame Fehler ist fast fünfmal höher als der nominell gewählte. Die Wahrscheinlichkeit, dass mindestens ein Vergleich zu einem stat. Fehlers 1. Art führt, beträgt 23%.
Lösungen zur Kontrolle der FWER
Paarweise Tests
- Bonferroni-Korrektur
- Kann immer angewendet werden, jedoch recht konservativ → Verlust an Teststärke
- Leicht zu berechnen: \(\alpha_{adj} = \frac{\text{Signifikanzniveau}~\alpha}{\text{Anzahl der durchgeführten Tests}} \Rightarrow \frac{0.05}{5} = 0.01\)
- Alternative Variante: Bonferroni-Holm-Verfahren (Sequentielle Bonferroni-Methode mit höherer Teststärke)
- Funktion
p.adjust()wandelt unkorrigierte p-Werte in korrigierte mittels Bonferroni-, Bonferroni-Holm- und weiteren Korrekturen um.
- Spezifische Tests für multiple Mittelwertvergleiche (sog. Post-hoc-Tests):
- Tukey’s, SNK, Ryan’s etc.
- passen Signifikanzniveaus auf unterschiedliche Weise an
Gruppenweise: Mehrstichprobentest
- Parametrisch:
- Varianzanalyse (ANOVA)
- Nicht-parametrisch:
- Kruskal-Wallis-Test (H-Test)
- Friedman-Test (Friedman-Rangsummentest)
Verteilungsfreie Tests
Kruskual Wallis Test (H Test)
- Logische Erweiterung des Mann Whitney U-Tests, Messung kann also auch auf einer Ordinalskala erfolgen
- Rang-basiertes Gegenstück zur 1-faktoriellen ANOVA → allerdings nur 95% der Teststärke
- Stichprobengröße der Gruppen kann unterschiedlich sein
- Annahme: Stichproben sollten annähernd die gleiche Verteilung haben
Friedman Test
- Rang-basiertes Gegenstück zur 1-faktoriellen ANOVA mit wiederholten Messungen
- Messung kann auch auf einer Ordinalskala erfolgen
- Anwendbar auf vollständige Blockdesigns
- Teststärke oder Effizienz: abhängig vom k/n-Verhältnis
Teststatistik H
\[H=\frac{12}{N*(N+1)}*\left( \sum\frac{R_i^2}{n_i} \right)-3*(N+1)\]
N = Gesamter Stichprobenumfang, \(n_i\) = Anzahl Messung pro Gruppe \(i\), \(R_i\) = Summe der Ränge pro Gruppe \(i\)
Teststatistik \(F_r\)
\[F_r=\left( \frac{12}{n*k*(k+1)}*\left( \sum R_i^2\right)\right)-3*n*(k+1)\]
n = Anzahl der Blöcke, k = Anzahl der Gruppen, \(R_i\) = Summe der Ränge pro Gruppe \(i\)
Kruskual Wallis Test
Einzelschritte
- Die Rohdaten aus den k Stichproben/Gruppen werden zu N Elementen zusammengefasst und dann vom niedrigsten zum höchsten Rang geordnet.
- Gleiche Werte werden mit dem Mittelwert der jeweiligen Ränge versehen.
- Diese Ränge werden dann wieder in die k Gruppen einsortiert.
- Bildung der Rangsummen pro Stichprobe i=1 bis k
- Berechnung der Teststatistik H
- Wenn mehr als 25 % den gleichen Rang haben, muss die H-Statistik korrigiert werden
- Prüfverteilung
- wenn \(n_i\) < 6 und k = 3 wird die H-Test-Tabelle verwendet
- wenn \(n_i\) ≥ 6 und k > 3 ist H näherungsweise \(\chi^2\)-verteilt (df = k-1)
Kruskal-Wallis-Test in R: manuell
Beispiel iris
Rangsummenbildung
df <- iris %>% mutate(ranks_sl = rank(Sepal.Length, ties.method = "average")) %>% group_by(Species) %>% summarise(sum_ranks = sum(ranks_sl)) df
# A tibble: 3 x 2 Species sum_ranks <fct> <dbl> 1 setosa 1482 2 versicolor 4132. 3 virginica 5710.
Teststatistik H und p-Wert
(h_sl <- 12/(150*(150+1))*sum(df$sum_ranks^2/50) - 3*(150+1))
[1] 96.76096
pchisq(h_sl, df = 3-1, lower.tail = FALSE)
[1] 9.741475e-22
Kruskal-Wallis-Test in R: kruskal.test()
Beispiel iris
Vektorsyntax
# x = abhängige Variable, g = Gruppierungsvariable/Faktor kruskal.test(x = iris$Sepal.Length, g = iris$Species)
Formelsyntax
# kruskal.test(x ~ g, data) kruskal.test(Sepal.Length ~ Species, data = iris)
Kruskal-Wallis rank sum test
data: Sepal.Length by Species
Kruskal-Wallis chi-squared = 96.937, df = 2, p-value < 2.2e-16Friedman Test in R
friedman.test()
Syntax
friedman.test(y, groups, blocks) # y = abhängige Variable # groups = Messwiederholung -> eigentlich zu untersuchender Faktor # blocks = zu untersuchende Objekte # oder Formelsyntax friedman.test(y ~ groups|blocks, data)
Your turn …
Quiz 1 zum Nachmachen
Kronblattweite in iris
Es soll die Hypothese getestet werden, dass sich die durchschnittliche Kronblattweite der drei Iris Arten unterscheidet.
kruskal.test(Petal.Width ~ Species, data = iris)
Kruskal-Wallis rank sum test
data: Petal.Width by Species
Kruskal-Wallis chi-squared = 131.19, df = 2, p-value < 2.2e-16
Varianzanalyse (ANalysis Of VAriance)
Was ist die Varianzanalyse oder ANOVA?
- In seiner einfachsten Form die Erweiterung der t-Statistik für mehrere Stichproben.
- Eine weit verbreitete statistische Methode.
- Testet auf Mittelwertunterschiede, indem die geschätzte Gesamtvariabilität in den Daten in 2 Teile zerlegt wird:
- Variabilität innerhalb (‘within’)der verschiedenen Gruppen
- Variabilität zwischen (‘between’) den verschiedenen Gruppen
- Wenn sich die Varianzen zwischen den Stichproben unterscheiden, müssen sich auch die Mittelwerte zwischen den Stichproben unterscheiden!
- Mit anderen Worten: je größer der Unterschied in der Varianz zwischen den Stichproben im Vergleich zum Unterschied innerhalb der Stichproben ist, desto wahrscheinlicher ist ein Unterschied im Mittelwert zwischen den Stichproben!
- Wichtige Annahme: Homogenität der Varianzen innerhalb jeder Gruppe (wichtiger als Normalverteilung!!!)
Variablen - und ANOVA-Typen
- Variablentypen:
- Unabhängige Variablen (X) sind meist kategorial und werden “Faktoren” genannt, die verschiedenen Kategorien werden “Faktorstufen” oder “Gruppen” genannt → sollten auch in R explizit als
factordefiniert werden! - Abhängige Variablen (Y) sind diskret oder kontinuierlich
- Unabhängige Variablen (X) sind meist kategorial und werden “Faktoren” genannt, die verschiedenen Kategorien werden “Faktorstufen” oder “Gruppen” genannt → sollten auch in R explizit als
- ANOVA-Typen (abhängig vom experimentellen Design):
- 1 Faktor = einfaktorielle Varianzanalyse (1-way ANOVA)
- 2 Faktoren = zweifaktorielle Varianzanalyse (2-way ANOVA)
- ≥3 Faktoren = mehrfaktorielle Varianzanalyse (multi-way ANOVA)
- Wenn >1 Faktor
- geschachteltes (nested) vs. gekreuztes (crossed) Design
- fixed vs. random vs. mixed effect model
- mit/ohne Messwiederholung, Blockdesign
Das Regressionsmodell
Das ANOVA Modell
Das Prinzip der ANOVA
1
Das Prinzip der ANOVA
1
Varianzzerlegung
Ein Anwendungsbeispiel: Lachszucht
Gewichtsunterschiede von Atlantischem Lachs, der in 4 verschiedenen Typen von Netzkäfigen gezüchtet wurde.
- Einfluss der Käfigtypen auf das mittlere Gewicht des Atlantischen Lachses (Salmo salar)
- Antwortvariable Y: Körpergewicht [in kg]
- Prädiktorvariable X: Käfigtyp
- kategorial mit 4 Gruppen: A,B,C,D
- Replikate sind Individuen pro Käfig i (\(n_i\) = 24)
Erkärte vs. nicht erklärte Varianz
| Replikat | Typ A | Typ B | Typ C | Typ D |
|---|---|---|---|---|
| 1 | 2.0 | 2.7 | 2.9 | 2.2 |
| 2 | 2.7 | 2.0 | 1.9 | 2.0 |
| 3 | 2.2 | 3.8 | 2.7 | 3.1 |
| … | … | … | … | … |
| 24 | 3.1 | 1.5 | 2.5 | 2.7 |
| Mittelwert | 2.5 | 2.5 | 2.5 | 2.5 |
Gesamtmittelwert: 2.5
- Behandlungen (hier = Käfigtypen) erklären nichts
- SSGruppen = 0
Erkärte vs. nicht erklärte Varianz
| Replikat | Typ A | Typ B | Typ C | Typ D |
|---|---|---|---|---|
| 1 | 5.0 | 4.1 | 2.9 | 0.8 |
| 2 | 5.0 | 4.1 | 2.9 | 0.8 |
| 3 | 5.0 | 4.1 | 2.9 | 0.8 |
| … | … | … | … | … |
| 24 | 5.0 | 4.1 | 2.9 | 0.8 |
| Mittelwert | 5.0 | 4.1 | 2.9 | 0.8 |
Gesamtmittelwert: 2.5
- Behandlungen (hier = Käfigtypen) erklären alles
- SSResidual = 0
Typische ANOVA Tabelle
| Sources of variation | SS | df | MS | F | p |
|---|---|---|---|---|---|
| Groups | \(n \sum(\bar{y_i}-\bar{y})^2\) | \(p-1\) | \(\frac{SS_{Groups}}{(p-1)}\) | \(\frac{MS_{Groups}}{MS_{Residuals}}\) | |
| Residuals (or within) | \(\sum\sum(y_{ij}-\bar{y_i})^2\) | \(p(n-1)\) | \(\frac{SS_{Groups}}{(p-1)}\) | ||
| Groups | \(n \sum(y_{ij}-\bar{y})^2\) | \(pn-1\) |
H0:
- Kein Unterschied zwischen den Gruppenmitteln \(\mu_1 = \mu_2 = ... = \mu_i = \mu\)
- Denn Gruppeneffekte sind definiert als \(\alpha_i = \mu_i - \mu\),
- Kann auch ausgedrückt werden als: “Es gibt keinen Effekt bestimmter Behandlungsstufen (treatment effects), die Unterschiede sind nur zufällig”: \(\alpha_1 = \alpha_2 = \alpha_3 = .. = 0\)
Sowohl \(MS_{Groups}\) als auch \(MS_{Residuals}\) schätzen \(\sigma^2\) → Quotient = 1
- F-Ratio folgt F-Verteilung mit \(df_{Groups}\) und \(df_{Residuals}\)
Beispiel: Gewichtsunterschiede bei Atlantischem Lachs
| Sources of variation | SS | df | MS | F | p |
|---|---|---|---|---|---|
| Groups | 52.10 | 4-1 | 171724 | 17.367 | <2e-16 |
| Residuals (or within) | 26.02 | 4*(24-1) | 0.283 | ||
| Groups | 78.12 | 4*24-1 |
- Es gibt einen signifikanten Effekt von mindestens einem Käfigtyp auf das mittlere Gewicht des Atlantischen Lachses. Mit anderen Worten: Die mittleren Gewichte von Atlantischem Lachs, der in verschiedenen Käfigen gezüchtet wurde, unterscheiden sich mindestens bei einem Käfigtyp.
- ABER: wir wissen noch nicht, welche Gruppe sich unterscheidet!
ANOVA in R
aov()
aov(formula, data = NULL, projections = FALSE, qr = TRUE, contrasts = NULL, ...)
Die volle ANOVA Tabelle erhältst Du aber nur mit der summary() Funktion, die auf das aov Objekt angewendet wird:
model <- aov(formula, data) summary(model)
ANOVA in R
Lachsbeispiel
mod <- aov(formula = weight ~ cages, data = salmon) mod
Call:
aov(formula = weight ~ cages, data = salmon)
Terms:
cages Residuals
Sum of Squares 53.91378 18.51342
Deg. of Freedom 3 92
Residual standard error: 0.4485898
Estimated effects may be unbalanced
summary(mod)
Df Sum Sq Mean Sq F value Pr(>F)
cages 3 53.91 17.971 89.31 <2e-16 ***
[ erreichte getOption("max.print") -- Zeile 1 ausgelassen ]
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Überprüfung der Modellannahme posteriori
op <- par(mfrow = c(2,2)) plot(mod)
par(op)
Effektgröße \(\eta^2\) bei der 1-faktoriellen ANOVA
plot.design()
Die Berechnung der Effektgröße für Designs zwischen Gruppen ist viel einfacher als für Designs innerhalb von Gruppen:
\[\eta^2 = \frac{SS_{Groups}}{SS_{Total}}\]
mod_sum <- summary(mod) # str(mod_sum) eta_sq <- mod_sum[[1]]$`Sum Sq`[1] / (mod_sum[[1]]$`Sum Sq`[1] + mod_sum[[1]]$`Sum Sq`[2] ) eta_sq
[1] 0.7443858
Effektgröße grafisch darstellen
plot.design()
plot.design(weight ~ cages, data = salmon)
Fragen?
Bei weiteren Fragen: saskia.otto(at)uni-hamburg.de
This work is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License except for the borrowed and mentioned with proper source: statements.
Image on title and end slide: Section of an infrared satellite image showing the Larsen C ice shelf on the Antarctic Peninsula - USGS/NASA Landsat: A Crack of Light in the Polar Dark, Landsat 8 - TIRS, June 17, 2017 (under CC0 license)